********************************************************************************
* Household level analysis: Income distributions
********************************************************************************

sum inc_lab_hh [fw = asecwth], d
sum inc_tot_hh [fw = asecwth], d

* Labor income distribution
_pctile inc_lab_hh [fw = asecwth], p(20,40,60,80,90)
scalar p20_inc_lab_hh = r(r1)
scalar p40_inc_lab_hh = r(r2)
scalar p60_inc_lab_hh = r(r3)
scalar p80_inc_lab_hh = r(r4)
scalar p90_inc_lab_hh = r(r5)
sum(inc_lab_hh) if inc_lab_hh<=p20_inc_lab_hh [fw = asecwth]
scalar mean_0_20_inc_lab_hh = r(mean)
sum(inc_lab_hh) if inc_lab_hh>p20_inc_lab_hh & inc_lab_hh<=p40_inc_lab_hh [fw = asecwth]
scalar mean_20_40_inc_lab_hh = r(mean)
sum(inc_lab_hh) if inc_lab_hh>p40_inc_lab_hh & inc_lab_hh<=p60_inc_lab_hh [fw = asecwth]
scalar mean_40_60_inc_lab_hh = r(mean)
sum(inc_lab_hh) if inc_lab_hh>p60_inc_lab_hh & inc_lab_hh<=p80_inc_lab_hh [fw = asecwth]
scalar mean_60_80_inc_lab_hh = r(mean)
sum(inc_lab_hh) if inc_lab_hh>p80_inc_lab_hh [fw = asecwth]
scalar mean_80_100_inc_lab_hh = r(mean)
sum(inc_lab_hh) if inc_lab_hh>p90_inc_lab_hh [fw = asecwth]
scalar mean_90_100_inc_lab_hh = r(mean)
sum(inc_lab_hh) [fw = asecwth]
scalar mean_inc_lab_hh = r(mean)

xtile inc_lab_hh_quintiles = inc_lab_hh [fw = asecwth], nquantiles(5)
tab inc_lab_hh_quintiles [fw = asecwth], matcell(x)
scalar pop_share_0_20 = x[1,1]/(x[1,1]+x[2,1]+x[3,1]+x[4,1]+x[5,1])
scalar pop_share_20_40 = x[2,1]/(x[1,1]+x[2,1]+x[3,1]+x[4,1]+x[5,1])
scalar pop_share_40_60 = x[3,1]/(x[1,1]+x[2,1]+x[3,1]+x[4,1]+x[5,1])
scalar pop_share_60_80 = x[4,1]/(x[1,1]+x[2,1]+x[3,1]+x[4,1]+x[5,1])
scalar pop_share_80_100 = x[5,1]/(x[1,1]+x[2,1]+x[3,1]+x[4,1]+x[5,1])
clear matrix

xtile inc_lab_hh_deciles = inc_lab_hh [fw = asecwth], nquantiles(10)
tab inc_lab_hh_deciles [fw = asecwth], matcell(x)
scalar pop_share_90_100 = x[10,1]/(x[1,1]+x[2,1]+x[3,1]+x[4,1]+x[5,1]+x[6,1]+x[7,1]+x[8,1]+x[9,1]+x[10,1])
clear matrix

scalar share_0_20_inc_lab_hh = pop_share_0_20*mean_0_20_inc_lab_hh/mean_inc_lab_hh
scalar share_20_40_inc_lab_hh = pop_share_20_40*mean_20_40_inc_lab_hh/mean_inc_lab_hh
scalar share_40_60_inc_lab_hh = pop_share_40_60*mean_40_60_inc_lab_hh/mean_inc_lab_hh
scalar share_60_80_inc_lab_hh = pop_share_60_80*mean_60_80_inc_lab_hh/mean_inc_lab_hh
scalar share_80_100_inc_lab_hh = pop_share_80_100*mean_80_100_inc_lab_hh/mean_inc_lab_hh
scalar share_90_100_inc_lab_hh = pop_share_90_100*mean_90_100_inc_lab_hh/mean_inc_lab_hh

di "Labor income distribution"
di share_0_20_inc_lab_hh
di share_20_40_inc_lab_hh
di share_40_60_inc_lab_hh
di share_60_80_inc_lab_hh
di share_80_100_inc_lab_hh
di share_90_100_inc_lab_hh

scalar drop _all

* Total income distribution
_pctile inc_tot_hh [fw = asecwth], p(20,40,60,80,90)
scalar p20_inc_tot_hh = r(r1)
scalar p40_inc_tot_hh = r(r2)
scalar p60_inc_tot_hh = r(r3)
scalar p80_inc_tot_hh = r(r4)
scalar p90_inc_tot_hh = r(r5)
sum(inc_tot_hh) if inc_tot_hh<=p20_inc_tot_hh [fw = asecwth]
scalar mean_0_20_inc_tot_hh = r(mean)
sum(inc_tot_hh) if inc_tot_hh>p20_inc_tot_hh & inc_tot_hh<=p40_inc_tot_hh [fw = asecwth]
scalar mean_20_40_inc_tot_hh = r(mean)
sum(inc_tot_hh) if inc_tot_hh>p40_inc_tot_hh & inc_tot_hh<=p60_inc_tot_hh [fw = asecwth]
scalar mean_40_60_inc_tot_hh = r(mean)
sum(inc_tot_hh) if inc_tot_hh>p60_inc_tot_hh & inc_tot_hh<=p80_inc_tot_hh [fw = asecwth]
scalar mean_60_80_inc_tot_hh = r(mean)
sum(inc_tot_hh) if inc_tot_hh>p80_inc_tot_hh [fw = asecwth]
scalar mean_80_100_inc_tot_hh = r(mean)
sum(inc_tot_hh) if inc_tot_hh>p90_inc_tot_hh [fw = asecwth]
scalar mean_90_100_inc_tot_hh = r(mean)
sum(inc_tot_hh) [fw = asecwth]
scalar mean_inc_tot_hh = r(mean)

xtile inc_tot_hh_quintiles = inc_tot_hh [fw = asecwth], nquantiles(5)
tab inc_tot_hh_quintiles [fw = asecwth], matcell(x)
scalar pop_share_0_20 = x[1,1]/(x[1,1]+x[2,1]+x[3,1]+x[4,1]+x[5,1])
scalar pop_share_20_40 = x[2,1]/(x[1,1]+x[2,1]+x[3,1]+x[4,1]+x[5,1])
scalar pop_share_40_60 = x[3,1]/(x[1,1]+x[2,1]+x[3,1]+x[4,1]+x[5,1])
scalar pop_share_60_80 = x[4,1]/(x[1,1]+x[2,1]+x[3,1]+x[4,1]+x[5,1])
scalar pop_share_80_100 = x[5,1]/(x[1,1]+x[2,1]+x[3,1]+x[4,1]+x[5,1])
clear matrix

xtile inc_tot_hh_deciles = inc_tot_hh [fw = asecwth], nquantiles(10)
tab inc_tot_hh_deciles [fw = asecwth], matcell(x)
scalar pop_share_90_100 = x[10,1]/(x[1,1]+x[2,1]+x[3,1]+x[4,1]+x[5,1]+x[6,1]+x[7,1]+x[8,1]+x[9,1]+x[10,1])
clear matrix

scalar share_0_20_inc_tot_hh = pop_share_0_20*mean_0_20_inc_tot_hh/mean_inc_tot_hh
scalar share_20_40_inc_tot_hh = pop_share_20_40*mean_20_40_inc_tot_hh/mean_inc_tot_hh
scalar share_40_60_inc_tot_hh = pop_share_40_60*mean_40_60_inc_tot_hh/mean_inc_tot_hh
scalar share_60_80_inc_tot_hh = pop_share_60_80*mean_60_80_inc_tot_hh/mean_inc_tot_hh
scalar share_80_100_inc_tot_hh = pop_share_80_100*mean_80_100_inc_tot_hh/mean_inc_tot_hh
scalar share_90_100_inc_tot_hh = pop_share_90_100*mean_90_100_inc_tot_hh/mean_inc_tot_hh

di "Total income distribution"
di share_0_20_inc_tot_hh
di share_20_40_inc_tot_hh
di share_40_60_inc_tot_hh
di share_60_80_inc_tot_hh
di share_80_100_inc_tot_hh
di share_90_100_inc_tot_hh

scalar drop _all